Mon. Not. R. Astron. Soc. OOP. [TlfTSl C2002') Printed 15 May 2009 (MN WI^ style file v2.2) 



o 
o 



Kinematics and history of the solar neighbourhood 
revisited 



Michael Aumer^* and James J. Binney^f 

Universitdtssternwarte Miinchen, Scheinerstr. 1, 81679 Miinchen, D 
^ Rudolf Peierls Centre for Theoretical Physics, University of Oxford, Oxford 0X1 3NP, UK 



Accepted . Received ; in original form 



in 
< 

6 



> 

(N 

in 

(N 

in 
o 
a^ 
o 



ABSTRACT 

We use proper motions and parallaxes from the new reduction of Hipparcos data and 
Geneva- Copenhagen radial velocities for a complete sample of ~ 15 000 main-sequence 
and subgiant stars, and new Padova isochrones to constrain the kinematics and star- 
formation history of the solar neighbourhood. We rederive the solar motion and the 
structure of the local velocity ellipsoids. When the principal velocity dispersions are 
assumed to increase with time as , the index (3 is larger for aw {Pw ~ 0.45) than 
for au {Pu « 0.31). For the three-dimensional velocity dispersion we obtain P = 0.35. 
We ex clude saturation of disc heating after ~ 3 Gyr as proposed bv lQuillen fc GarnettI 
(|2000t ). Saturation after > 4 Gyr combined with an abrupt increase in velocity disper- 
sion for the oldest stars cannot be excluded. For all our models the star-formation rate 
is declining, being a factor 2-7 lower now than it was at the beginning. Models in which 
the SFR declines exponentially favour very high disc ages between 11.5 and 13 Gyr 
and exclude ages below ~ 10.5 Gyr as they yield worse fits to the number density and 
velocity dispersion of red stars. Models in which the SFR is the sum of two declin- 
ing exponentials representing the thin and thick discs favour ages between 10.5 and 
12 Gyr with a lower limit of ^ 10.0 Gyr. Although in our models the star-formation 
rate peaked surprisingly early, the mean formation time of solar- neighbourhood stars is 
later than in ab-initio models of galaxy formation, probably on account of weaknesses 
in such models. 
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1 INTRODUCTION 

For the solar neighbourhood we have the most detailed ob- 
servational data available for any galactic disc. In particular, 
the distributions within the solar neighbourhood of stellar 
ages, metallicities and space velocities are keys to deduc- 
ing how the disc has evolved chemically and dynamically. 
Work directed at understanding how galaxies have formed 
and evolved is a major area of contemporary astronomy, and 
studies of the local disc have an important role to play in 
this effort. Despite much progress many questions about the 
evolution of the local disc remain open. 

The Hipparcos catalogue provides uniquely useful data 
for studies of the solar neighbourh ood, because the data are 
homogeneous and of high quality. iDehnen fc BinnevI l| 19981 . 
hereafter DB98) defined a kinematically unbiased sample 
of Hipparcos stars and investigated its kinematics, while 
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iBinnev. Dehnen fc BerteiiH (|2000l . hereafter BDBOO) used 
the sample to model the history of star formation and the 
stochastic acceleration of stars in the disc. In this paper we 
redefine the sample and then re-work these papers for sev- 
eral reasons. 

(i) In a systematic re-redu ction of the astrom etric data 
from the Hipparcos satellite, Ivan LeeuwenI (|2007l ) has been 
able to diminish significantly the errors for a large number 
of stars. Using the new data we can enlarge the sample of 
stars that can be used for modelling. 

(ii) BDBOO did not take into account the variation of the 
scale height of stars with age. Consequently, their con- 
clusions regarding past star-formation rates and rates of 
stochastic acceleration are wrong. Our analysis remedies this 
defect. 

(iii) From an analysis of the ag e-velocity dispersion rela - 
tion in the sample o f 189 stars in lEdvardsson et akl (|l993l ) , 
iQuillen fc GarnettI (|2000l ) argued that stochastic heating 
saturates after only 3 Gyr and that the velocity dispersion 
of disc stars increased abruptly ~ 9 Gyr ago as a result of a 
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Figure 1. Colour- magnitude diagram of the magnitude limited 
subsample of 20 360 Hipparcos stars with relative parallax er- 
rors smaller than 10%. The lines indicate our selection of main- 
sequence stars, there are 15113 stars between the lines 



merger. We use our much larger sample to test this conjec- 
ture. 

(iv) Metallicity and radial-v elocity measurements fr om the 
Geneva Cop enhagen Survey iNordstrom et al.l l|2004l ). here- 
after GCS, [H olmberg et al.1 (|2007l ). hereafter GCS2 and 
iHolmberg et al. (2008.), hereafter GCS3] of F and G dwarfs 
allow more accurate determinations of the metallicity dis- 
tribution of disc stars and the three-dimensional velocity 
dispersions in the corresponding colour interval. 

(v) New Padova is ochrones have recently been published 
l|Bertelli et al. I [20081 ). allowing us to update the input con- 
cerning stellar evolution and to improve the modelling of the 
metallicity distribution. 

(vi) A revision of the impact of interstellar reddening on our 
sample changes some results significantly. 



The goal of this paper is the same as that of BDBOO: 
to understand as much of the history of the disc as we can 
by modelling the large-scale structure of velocity space near 
the Sun. We do not address the small-scale structure of this 
space, which was first clearl y revealed by the H ippar cos Cat- 
alogue through the work of ICreze et al.i (|l998l ) and iDehnenI 
lll998t) . 

Section [5] describes our input data. Section |3] updates 
DB98 by extracting the kinematics of the solar neighbour- 
hood. The work of BDBOO is updated in Section [l] which 
describes our models, and Section [S] which describes fits of 
them to the data. Section [6] sums up and compares our re- 
sults with those of other authors. 



2 THE DATA 
2.1 The Sample 

We follow the procedure of DB98 to select from the new 
reduction of the Hipparcos catalogue (|van LeeuwenI [2OO7I ) 
a kinematically unbiased, magnitude-limited subsample of 
single stars with high-quality astrometric data. 

We determine the apparent magnitude up to which 
the Hipparcos catalogue is complete (Viim ^ 8 mag), which 
depends on colour and position on the sky, by compar- 
ing the Hipparc os catalogue with the Tycho2 catalogue 
l|H0g et al.ll200d ). which is complete to about Vrycho ~ H- 
For 16 X 16 X 10 uniformly spaced bins in galactic coordinates 
sin 6 and I and colour {B — F)Tycho in (—0.3, 1.5), we select 
those stars that are brighter than the second brightest star 
per bin that is in the Tycho2 but not the Hipparcos cata- 
logue. 44 567 out of the 118 218 Hipparcos stars are single 
stars that fulfil this criterion. From these stars we further 
select stars with relative parallax errors of 10% or less. This 
criterion reduces the sample from 44 567 to 20 360 but en- 
sures that the measured proper motions yield fairly accurate 
tangential velocities. 

Fig.[T] is the colour-magnitude diagram of the sample. 
For our analysis of the star-formation history, we require a 
one-to-one relationship between the maximum lifetime of a 
star and its colour, so we restrict ourselves to main-sequence 
stars. The lines within which stars are deemed to lie on the 
main sequence are shown in Fig.[T] The applied CMD cut 
leaves the possibility for a small number of subgiant stars to 
enter the sample in the colour range 0.6<i3 — y<0.75, 
where the sample might thus be slightly biased towards old 
stars. The 15 113 stars that lie between these lines comprise 
our final sample, which is approximately 27% larger than 
that of DB98. The growth in the sample is most pronounced 
for blue stars since these tend to be more luminous and 
distant and therefore have the smallest parallaxes. 

GCS have measured radial velocities for 6 918 single F 
and G dwarfs in our sample. These stars are confined to the 
colour interval 0.4 < B — V < 0.8. We use this subsample, 
for which individual space velocities can be determined, as 
a control of the results obtained with the main sample. 



2.2 Isochrones and Metallicities 

The age distribution of main-sequence stars of a given colour 
must vary with colour because at the blue end of the main 
sequence all stars must be younger than the short main- 
sequence lifetime there, while at the red end of the main 
sequence we see stars that are as old as the disc in ad- 
dition to recently-fo rmed stars. We use Padova isochrones 
(|Bertelh et al.ll2008l ) for masses ranging from 0.15 to 4M0 
to determine the age distribution at each colour. Isochrones 
depend significantly on the metal content Z and the Helium 
abundance F, so we have to use several isochrones to simu- 
late the metallicity distribution of the solar neighbourhood. 
We refer to Section [5TT] for a discussion of the weightings of 
the isochrones employed and focus on which isochrones to 
employ. 

For the conversion from solar u nits to Z and Y we 
need to know the solar abundances. iBertelh et al] (|2008l ) 
argue that the update of the solar chemical composition 
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Table 1. The metallicities of the Padova isochrones that were 
used in the models and the corresponding [Fc/H] values, which 
were calculated with solar abundances Zq = 0.017 and Yq = 
0.260. 



z 


Y 


[Fo/H] 


0.002 


0.23 


-0.95 


0.003 


0.23 


-0.78 


0.004 


0.23 


-0.65 


0.006 


0.24 


-0.47 


0.008 


0.24 


-0.34 


0.010 


0.25 


-0.24 


0.012 


0.25 


-0.16 


0.014 


0.25 


-0.09 


0.017 


0.26 


0.00 


0.020 


0.27 


0.08 


0.026 


0.28 


0.20 


0.036 


0.30 


0.36 



by iGrevesse et al] (|200'if ) (Yq = 0.2486,^0 = 0.0122) 
gives rise to several problems and uncertainties and adopt 
(■^0 = 0.260, Z0 = 0.01 7) for their solar model. More- 
over IChaplin et all (|2007h report that solar models with 
Z < 0.018 7 are inconsist e nt wi th the helioseismic data. Since 
we use the lBertelli et all (|2008l ) isochrones, we use their solar 
composition. 

As the measurement of the Helium abundance y of a 
star is difficult, there is little data available on the relation- 
ship Z{Y) within the solar neighbourhood. A common, but 
not undisputed procedure is to assume a linear enrichment 
law Y = Yp + {dY/dZ)Z , where Yp is the primordial he- 
lium abundance and AY/AZ is the helium-to- metal enrich- 
ment ratio. Jimenez et al. (,200a) find dY/AZ = 2.1 ± 0.4 for 



Yp = 0.236, whereas Casagrande et ahl (|2007l ) using Padova 



isochrones find that it is doubtful whether a universal lin- 
ear enrichment law exists. For metallicities similar to the 
solar one they also find AY/AZ = 2.1 ± 0.9. Unfortunately, 
this slope when combined with either the traditional values 
{ Yq,Zq) = (0.260 , 0.017 ) or the new values (0.2486,0.0122) 
ofE 

revesse et al.l |2003) is inconsiste nt with the generall y 
accepted WMAP value = 0.24815 (|Spergel et al.ll2007l ). 
Nonetheless we adopt AY/AZ = 2.1. Table[T]gives the chem- 
ical compositions of the isochrones we have employed. 



3 STELLAR KINEMATICS 

We follow the method of DB98 to calculate the three- 
dimensional velocity dispersion cr as a function of colour 
and use a sliding window in B — V with borders that were 
adjusted to ensure that there are always 500 stars in a win- 
dow. Every 100 stars a new point is plotted and thus every 
fifth point is statistically independent of its predecessors. 
As explained in DB98, the averages we use here are sen- 
sitive to outliers, which have to be rejected. We therefore 
use an iterative method and reject stars that contribute to 
= air + ay + o\r more than times the value of the 



previous iteration^ We use k = 3.5 and 6 iterations; 55 stars 
were rejected. 

The red points in Fig. [2] show the resulting run of a 
against the mean B — V colour in the sliding windows, with 
errors in o shown for every fifth point. Parenago's discon- 
tinuity - the abr upt end to the i ncrease in velocity disper- 
sion with B — V l|Parenagolll950l ) - is beautifully visible at 
B — V~ 0.61. We also calculated a by combining GCS radial 
velocities with Hipparcos proper motions, and in Fig. [2] the 
results are plotted in blue. Where the samples are compara- 
ble, there is excellent agreement between the two measure- 
ments of (T. At the upper and lower limits of the GCS sample 
there are discrepancies because the GCS stars were selected 
by spectral type, so the bluest GCS bins only contain F 
stars, while the bins from the Hipparcos sample contain A 
and F stars and consequently have a lower mean age. For 
the reddest bins the situation is the same for G and K stars, 
but this time counterintuitively the pure G star sample has 
a higher mean age and velocity dispersion - we explain this 
phenomenon in Section [5.81 

The diagonal elements of the velocity-dispersion tensor 
are also of interest, as the stochastic acceleration mecha- 
nisms at work are anisotropic so each component of the ten- 
sor can evolve independently in time (cf. e.g. GCS). Fig.|3] 
plots (Ju and ay against the mean B — V colour of each bin 
and Fig. 13] shows aw- The three-dimensional dispersion is 
dominated by the radial dispersion au, which consequently 
shows the same features as a. The other two components 
are different, however, and Parenago's discontinuity is not 
as beautifully visible in them. For cry, the dispersion in the 
direction of rotation, we notice that the slope changes at 
B — V~ 0.6, but it continues to increase until B — V ~ 0.75. 
Interestingly, in Fig. [4] the values aw from proper motions 
(red points) show a significant bump a.t B — V ^ 0.6. A 
similar bump is visible in DB98 but less strikingly so on ac- 
count of the larger bins used by DB98. Values of aw from 
GCS space velocities have smaller errors and show a less 
distinctive bump. The GCS error bars on individual com- 
ponents of velocity dispersion are significantly smaller than 
the corresponding error bars from the larger Hipparcos sam- 
ple because each Hipparcos star in effect constrains only two 
components of velocity dispersion. The error bars from the 
two measurements of aw overlap, so the results are not in- 
consistent. The bump originates from a decrease in average 
age redwards of the discontinuity that we will explain in 
Section [5.81 In subsequent work we require a functional fit 
to the data in Fig.|4l Since the radial velocities add valuable 
additional information, we fit a 5th-order polynomial to the 
GCS values where these are available, and elsewhere to the 
proper-motion values. The black curve in Fig. [4] shows the 
resulting relation. 

As in DB98, the afjy component of the velocity dis- 
persion tensor is nonzero, implying that the principal axes 
of a^ are not aligned with our coordinate axes. We diago- 
nalised a^ to find its eigenvalues al, the ratios of the square 
roots of which are plotted in Fig.[S] The ratio a\/ai (red 
and magenta) shows a slight decrease with increasing B — V 



^ We use a right-handed coordinate system with x increasing 
towards the Galactic centre and y increasing in the direction of 



Galactic rotation. 
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Figure 2. Total velocity dispersion vs. colour. The red line connects the data points for the van Leeuwen- ffipparcos proper motions and 
here and elsewhere only independent data points carry Icr-crrorbars. The blue line is the relation obtained from the GCS space velocities. 
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Figure 3. In-plane components of velocity dispersion vs. colour. 
The red line connects the data points for the ajj component and 
the black line connects the data points for the ay component, 
both calculated from the van Leeuwen- ffipparcos proper motions 
with bins of 1000 stars. We also show cry and ay from the GCS 
space velocities with bins of 500 stars. Only independent data 
points carry Icr-errorbars. 
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Figure 4. Vertical velocity dispersion vs. colour. The red line 
connects data points for aw from the van Leeuwen— ffipparcos 
proper motions with bins of 1000 stars, while the blue points 
show values from the GCS space velocities with bins of 500 stars. 
Only independent data points carry Itr-errorbars. The black line 
shows the polynomial fit that is explained in the text and used in 
subsequent work. 



and therefore age, dropping from ~ 1.9 to ~ 1.6 through the {Uo, Vb, Wo) as explained in DB98 and found: 
range 0.3 < B ~V < 0.55. The ratio a\/a3 drops from ~ 2.8 _ gg ^ q 33) i^nis"^ 

for the bluest stars to ~ 2 (judging from the GCS sample) ' _^ 

for stars redder than Parenago's discontinuity. ^0 = (5.25 ± 0.54) km s (1) 

We calculated the sun's velocity with respect to the LSR Wa — (7.07 ± 0.34) km s~ . 
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Figure 5. The ratios of the eigenvalues of the velocity dispersion 
tensor for bins with 1000 stars and a new bin after 200 stars. 
(Ti/(T2 from proper motions (red) and GCS space velocities (ma- 
genta). (Ti/(73 from proper motions (blue) and GCS space veloc- 
ities (black). 



These values are consistent with the results of DB98 but 
our er ror bars are smaller by ~ 11 per cent. Ivan LeeuwerJ 
((2003), using a sample of ~ 20 000 main-sequence stars se- 
lected by imposing only an upper limit for relative parallax 
errors of 10%, found Uo = 10.77 ± 0.36, Vo = 3.21 ± 0.52, 
in conflict with our values, and Wo = 7.04 ± 0.14, which 
is consistent with our findings. Since van Leeuwen's sample 
extends below the limit to which the Hipparcos sample is 
photometrically complete, his solar motion may be affected 
by the kinematic biases that are known to be present in 
the Hipparcos input catalogue: stars thought be interesting 
for a variety of reasons were added to the catalogue, and 
many such stars had come to the attention of astronomers 
by virtue of a high proper motion. 

Using our values for the solar motion and for a^^ to cal- 
culate the coefficient in Stromberg's asymmetric drift equa- 
tion 



(2) 



we find k = 74±5kms~^, which is lower than, but still con- 
sistent with the result of DB98, a nd also with the estimated 
scale length of the Galactic disc (|Binnev fc Tremainell200^ . 
§4.8.2(a)). 

As in DB98, we also calculated the vertex deviation 1^, 
the angle by which one has to rotate the applied Cartesian 
coordinate system around its z axis to make the velocity 
dispersion tensor diagonal in the (v^cVy) plane. Fig.|5]shows 
the vertex deviation as a function of colour B — V for bins 
containing 1000 stars each. As in DB98 the vertex deviation 
has a minimum at the colour of Parenago's discontinuity. 
The displayed errors for the Hipparcos sample are slightly 
larger than in DB98 because we use smaller bins. 
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Figure 6. The vertex deviation ly as a function of colour B — V 
for bins containing 1000 stars each. 



4 MODELLING THE HISTORY OF STAR 
FORMATION AND HEATING 

4.1 Interstellar reddening 

For our models, it is important to apply a correction for 
interstellar reddening to the data because reddening shifts 
the location of Parenago's discontinuity, which strongly in- 
fluences the results of our models for the star formation his- 
tory. BDBOO corrected the data assuming a linear reddening 
law, E{B — V) — 0.53 mag kpc~^. This law is appropriate for 
distances of order kiloparsecs, but we live in a low-density 
bubble within the ISM, so the column density of the ISM 
and thus the amount o f reddening d epends on the direction 
of the line of sight (e.g. iFrischll 19951 ) and for small distances 
is typic ally smaller than the mean relation would imply. The 
code of lHakkila et al.l l| 19971 ) allows us to determine the red- 
dening cis a function of galactic coordinate s (l,b), but unfor- 
tunate ly cannot reproduce th e results of iLallememt et al] 
(|2003l ) or lVergelv et al.l (|l998l ). who find that the reddening 
within ~ 70 pc of the sun is essentially negligible. Recently 
GCS2 found that reddening is negligible within 40 pc and 
E{b-y) = 0.0048 mag between 40 and 7Opc0 

We decided t o deredden according to Fig. 4 of 
IVergelv et al.l (jl998l ) . For stars within 40 pc of the midplane 
and with a distance from the sun between d = 70 and 300 pc, 
they find E{h ~ y) ^ 0.35(d — 70pc)/kpc, so we have used 





0.47(d - 



70 pc)/lkpc 



for d < 70 pc 
otherwise. 



(3) 



Fig-Eldisplays the effect of dereddening on the {B—V, a) 
diagram by comparing data dereddened according to the old 
and new prescriptions. The binsize used here was 1000 stars 
and every 200 stars a new data point was added. We see that 
BDBOO considerably overestimated the effect of reddening. 



^ The convers ion factor between colour s ystems is E{B — V) ! 
1.35 £(6 - y) l|Persinger fc CasteladFiggoh . so GCS2 find E{B ■ 



V) ~ 0.0065 between 40 and 70 pc. 
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Figure 7. The effect of dereddening on ttio a vs. {B — V) diagram 
for bin sizes of 1000 stars. 

4.2 Volume completeness 

As the velocity dispersion vs. colour C = B ~ V diagram 
alone allows models with a rather large variety of parame- 
ters, we use the number of stars per colour interval, AN /AC, 
as an additional constraint. Our model of the dynamics will 
predict N^yiiC), the number of stars of a given colour in a 
vertical cylinder through the disc that has radius Rcyi and 
the sun on its axis. Two points have to be borne in mind 
when relating this to the observed number density AN /AC: 
a) the scale height of stars varies with velocity dispersion and 
thus colour, and b) the radius of the sphere within which the 
sample is complete varies with colour. 

Assuming the stellar distribution is in equilibrium and 
neglecting variations with galactic radius, the stellar number 
density as a function of vertical coordinate z, potential <3?(z) 
and vertical velocity dispersion aw(C) is 



v{z, C) = vo{C) exp 



(4) 



where the normalising factor i^o(C') is to be determined from 
the observed density d7V/dC Integrating u[z,C) through 
the cylinder we have 



A^cyl(C) = Uo{C)tvR%, 



Az exp 



(5) 



Integrating v{z, C) through the sphere within which the star 
count is complete at the given colour we have 



Az {R^{C) — z ) exp 



(6) 



= vo{C)I{C), 

where Rc{C) is the completeness radius and zo is the vertical 
distance of the sun from the galactic midplane. We now have 
that the normalising constant vq is 

AN/AC 



I{C) 



(7) 



and on substituting this into equation ^ we obtain the re- 
quired relation between A'^cyi and the star counts AN /AC. 
Rcyi can be set to any value larger than the largest value 



100 




100 120 



distance in pc 

Figure 8. Distance histogram for the colour bin with 0.459 < 
B — V < 0.476. The dashed line shows the expected relation for 
a volume complete sample. The final radius used was 71pc and is 
indicated by the dotted line. 



of R c(C). We set zo = 15 pc (jBinnev et all Il997l : Ijoshil 
l2007f ) and for aw{C) we use the polynomial fit shown in 
Fig. in For the pote i itial ^ (z) we use Model 1 of §2.7 of 
iBinnev fc Tremaind (|2008D . 

For the determination of Rc we compare the distribu- 
tion of distances of stars within a radius 7? in a given colour 
bin to the distribution we would expect, namely 



n{R, C) = iyoiC)Tv / Az {R^ - z^) exp 



(8) 



The model density n{R,C) depends on Rc through i^o(C) 
and we adjust Rc until we have a value 0.4 < P < 0.6 
for the Kolmogorov-Smirnov (KS) probability that the ob- 
served distance distribution at R < Rc is consistent with 
the model distribution. These comparisons were made for 
colour bins that contain 500 stars each. Fig. [8] shows a typi- 
cal histogram and the corresponding model distribution (IS}. 
We finally applied a 5th order polynomial fit to the relation 
Rc{C) (cf. Fig.[9ll and used these radii in eqs. ([5]) and (O to 
relate A^cyi to AN /AC. An average of 49 per cent of the stars 
in a colour bin lie within RdC). 

A priori, it is not clear, how to choose the exact lim- 
its for the KS probability. Fig.[S] might give rise to doubts 
concerning the completeness out to the chosen value of Rc. 
However, if we require higher KS probabilities resulting in 
lower completeness radii and smaller numbers of stars avail- 
able for analysis, the results do not change in any significant 
way. We estimated the errors in the number of stars in the 
column by varying Rc by 10 per cent around our preferred 
value and by varying aw by its Icr-errors. The resulting er- 
ror was added in quadrature with the Poisson error. This 
procedure yielded relative errors that varied from below 10 
to 25 per cent, but we imposed a lower limit of 15 per cent 
on the relative error. 

To relate all this to BDBOO, Fig. llOl shows our selection 
function (blue) and that of BDBOO (red). One fundamental 
difference between BDBOO and this work is that we con- 
sider a column stretching to oo, whereas BDBOO considered 
a sphere with a radius of 100 pc. In BDBOO stars within 
this volume more luminous than 5I/q were assumed to be 
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Figure 9. The blue data points show Rc as determined by the 
Kolmogorov-Smirnov test. The red Hne displays the applied 5th 
order polynomial fit. 




0.1 r- 



0.01 z- 



0.001 

0.4 0.6 
B-V 

Figure 10. The selection function as in BDBOO (red) and in our 
models (blue) 



complete, and for less luminous stars the radius of complete- 
ness was supposed to decrease as L^^^. Hence in Fig. llOl the 
BDBOO selection function is unity bluer than B — V ^ 0.35 
and then declines linearly. Since we consider a column that 
extends to infinity, the selection function is always a declin- 
ing function of luminosity and therefore colour. Hence the 
most significant difi^erence between the two selection func- 
tions is at the blue end where (a) the ability of luminous 
stars to enter the sample even when far from the plane 
makes them relatively more numerous in the sample, and 
(b ) the m ore accurate reduction of the raw Hipparcos data 
by Ivan lie euwcii (2007) has made stars with higher distances 
available. Even at fainter magnitudes the introduction of 
the scale-height correction has slightly increased the slope 
of the selection function, thus depressing the chances of a 
red star to enter the sample. Consequently, our models have 
to increase the predicted numbers of red stars in a cylinder 
relative to the predictions of the BDBOO models. 

In summary, roughly half the sample stars contribute to 
the values of A^'cyi that we model, but all of the sample stars 
contribute to the modelled values of the velocity dispersions. 



4.3 Distribution over age of stars at a given colour 

Following BDBOO the distribution of main-sequence stars in 
a certain volume over age and mass is given by: 



d"iV ^ r C(M) SFR(t) for r < r^MM), 
dA/dr 1 otherwise, 



(9) 



where Tmax(M) denotes the main-sequence lifetime of a star 
of initial mass M, £,{M) the initial mass function and SFR(r) 
the star formation rate. 

BDBOO used a Salpeter-hke power-law IMF ^(M) oc 
M" and an exponential SFR oc exp(7r). They found that 
the corresponding character ictic p arameters a an d 7 were 
strongly correlated (see also lHavwood et allll997h . A high, 
positive value of 7, i.e. a higher SFR in the past, relatively 
increases the number of red stars compared to blue stars. 
This effect can be cancelled by a relatively fiat IMF creating 
more blue stars. 

Because of this co rrelation, we decided to use the IMF 
oflKrouEaelaD (|l993l ): 



C(Af) « 



A4-2'2 
Af" 



if 0.08 Mq < M < 0.50 M© 
if 0.50 Mq < M < LOOM© 
if 1.00 Mq < M < 00. 



(10) 



The power- law parameter a for M > 1.0 Mq has the 
strongest infiuence on our results and is thus allowed to vary 
around a — —2.7, the value found bv lKroupa et al. I ((19931). 

With the IMF fixed within a small range, we tested 
several models for the star formation history: 



(i) A simple exponential SFR 
SFR(t) (X exp(7r) 

(ii) A SFR of the form 
SFR(t) (X A exp(Ar) + exp(7r). 



(11) 



(12) 



with A > 7, adding an additional amount of star formation 
in the early universe, 
(in) A SFR of the form 



SFR(r) (X 



r2 — r 



(13) 



[{t2 - r)2 + fo2]2 ' 

as proposed by I Just fc Jahreis3 (120071 ). 

(iv) A smooth SFR o verlaid with a factor varyin g with time 
according to Fig. 8 of iRocha-Pinto et af] ()200(]| ). 

We need the isochrones that were described in Section 
12.21 to determine the mass range that can be found in a 
colour interval at a given time. We cut each isochrone ofi' 
above the point where it is 1.8 mag more luminous than the 
ZAMS at the same colour. The isochrones provide infor- 
mation only at a limited number of times, so at a desired 
time one has to interpolate between the next older and next 
younger isochrones; for the details of this we refer to BDBOO. 

We can relate A'^cyi to the average distribution over age 
in a colour interval (C-, C+) by 



/d7V\ _ 1 exp(7r) 

\ d7/(C_,c+) " C+-C- (1 + a) 



.(14) 



Here the sum is over all mass ranges {M-j,M+j) that lie 
in the colour interval at age r. 



8 Michael Aumer & James J. Binney 



4.4 Age— metallicity relation 

As we use isochrones with a significantly higher number 
of metaUicities than in BDBOO, we are able to include a 
variation of metallicity with age. As guidelines we use the 
metallicity distribution published in GCS2 and the age- 
met alli city relation from the models of ISchonrich fc BinnevI 
l|2009al . hereafter SB09a), who are able to reproduce the find- 
ings of GCS2 and whose r esults are also co nsistent with the 
age-met allicity relation of iHavwoodl (|2008l ). 

We find that the GCS2 metallicity distribution has 
mean [Fe/H] = —0.12 and the dispersion 0.17 d ex. This re- 
sult is similar to that of lGirardi fc SalarisI (1200 ll ) , who found 
that the metallicity distribution for K giants is well repre- 
sented by a Gaussian with a mean of —0.12 and a dispersion 
of 0.18 dex. IHavwoodl (|200lh proposed that the distribution 
was centred on solar metallicity, which shows that the un- 
certainty is not negligible. Using a Gaussian representation 
however omits the metal-poor tail of the GCS2 distribution, 
which comprises approximately 4% of the total sample and 
spreads in [Fe/H] from —1.2 to —0.4. 

From the models of SB09a we extracted that stars with 
[Fe/H] < —0.7 form only in the first ~ 0.6 Gyr, stars with 
—0.7 < [Fe/H] < —0.4 show a high contribution from 
the time interval 0.5 — 1.5 Gyr, but also a younger com- 
ponent and stars with [Fe/H] > —0.4 started forming after 
~ 0.6 Gyr and still form today. 

We therefore construct a model that features the fol- 
lowing two components: 

(I) The 'thin disc' component 

This component comprises ^ 96% of the model stars and 
is represented by a Gaussian distribution with the above 
mentioned mean, so its stars belong to the isochrones with 
the ten highest metaUicities from Table[T] The intrinsic dis- 
tribution of metaUicities will be narrower than the observed 
one on account of observational errors (cf. Section 15. ip . 
GCS2 conclude that there is no significant change with age 
in the mean metallicity of solar-neighbourhood stars, so in 
our models we can consistently use a fit to the current dis- 
tribution at all times. The models for the star-formation 
history as described in Section [4.31 and the age Tmax refer to 
this component only. 

(II) The 'low-metallicity' component 

For this component we use only the isochrones with the 
four lowest metaUicities from Table[T] 

Stars represented by the isochrones with [Fe/H] = —0.95 
and —0.78 with contributions of ~ 0.5% and ~ 1.0% have 

ages T G (Tmax, Tmax -|- 0.6 Gyr). 

Stars represented by [Fe/H] = —0.65 with a total con- 
tribution of ~ 1.5%; two thirds of these have ages t £ 
(Tmax — 0.5 Gyr, Tmax + 0.1 Gyr) and one third has ages 
T G (Tmax - 5 Gyr, Tmax — 0.5 Gyr). 

Stars represented by [Fe/H] = —0.47 with a contribution 
of ^ 1.0% have ages t G (Tmax — 1 Gyr, Tmax). There is also 
a significant contribution from stars with this metallicity to 
the 'thin disc' component (cf. Section [5. ip . 

It would seem desirable to assign a velocity dispersion as 
high as the one of the 'thick disc' to this component, however 
this significantly diminishes the quality of the fits. We thus 
model both components with a single disc heating rate as 



described in the following Section and discuss this result in 
Section 15.71 



4.5 The disc heating rate 

We model the velocity dispersion cr^ of a group of stars with 
a known distribution in age by 

2 ^ dT(djV/dT)a^(T) 
/;"-''dT(diV/dT) 

A sim ple and often used model (e.g. iBinnev fc Tremaind 
|2008| . §8.4), that we will mainly consider here is 



T -I- Tl 

10 Gyr + Ti 



(16) 



where «io and ti characterise the velocity dispersion at 10 
Gyr and at birth and /3 describes the efficiency of stochas- 
tic acceleration. However, iQuiUen fc GarnettI l|2000l ) argued 
for saturation of disc heating after ~ 3 Gyr and an abrupt 
increase in velocity dispersion at ~ 9 Gyr, which they con- 
nected to the formation of the thick disc. Therefore we also 
consider the following model 



a(T) = 



Vs 
Vs 



T + Tl 



for 

for 
for 



T < Tl 

Tl < T < Ta 
T > Ta, 



(17) 



where Ti and Ta are the times for the occurrence of satura- 
tion and the abrupt increase, Vs and riVs are the saturation 
velocity dispersions for the thin and the thick disc and (3 
and Tl are as above. 

Our velocity dispersion data from Section [3] are not for 
an infinite cylinder, but for a volume- limited Hipparcos sam- 
ple. As a coeval population heats, it will spread in z and the 
fraction of its stars that contribute to the Hipparcos sphere 
will drop. So the contribution of this population to the mea- 
sured dispersion will be less than that of a younger popu- 
lation. We resolve this problem by introducing a weighting 
factor F{t) ^ 1 such that a population's contribution to the 
measured dispersion is proportional to 

dN_ 

d-T heating 



Fir) 'I . 

(IT cylinder 



(18) 



To estimate F{t) we consider the conservation of the num- 
ber of stars born in a certain time interval dT in the ap- 
proximation that we can neglect radial mixing. Then, as the 
population heats and spreads in z, its central density drops 
by a factor F{t), which is given by 



F(t) 



dz exp 



constant. 



(19) 



The contribution of a population to any sphere around the 
sun can now be obtained from equations ^ and (O. For 
our first fit of the model to the data we take F{C,t) = 1. 
The resulting function awiT) is used to determine F{C,t), 
a new fit is made and F is redetermined. This sequence of 
operations is rapidly convergent. 

The Levenberg-Marqu ardt non-linear least-squares al- 
gorithm l|Press et al.l Il986l ) is used to minimise the of 
the fits to the data for a{B - V) and dN/dC{B - V). The 
parameters adjusted are (for the standard disc heating and 
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SFR model) a, (3, 7, rmax, n and v\_o. Data a.t B — V < Q 
are discarded for fear that the sample of young stars is kine- 
matically biased. As stated above, the stars are binned in 
sliding windows of 500 stars, a new one every 100 stars. It 
might seem desirable to use only statistically independent 
bins (every fifth bin), but then the results turn out to de- 
pend on which subset of bins is used. The compromise used 
was to use every third bin, which reduces the degradation of 
the information available about the colours at which disper- 
sions change, at the price of yielding values of that are 
slightly too low. 



5 THE RESULTS OF THE FITS 

As described in BDBOO, there are correlations between the 
parameters which limit the usefulness of the formal errors 
on the parameters. So we present the results in the follow- 
ing form: we first show the infiuence of different metallicity 
weightings on the results and after standardising on a plau- 
sible configuration demonstrate the possible range of each 
parameter by showing the results obtained for fixing it to 
certain values and leaving the other parameters free. We 
generally use the total velocity dispersion for our models. 
The components of a are studied in Section [5.51 

5.1 The influence of the metallicity weighting 

It is interesting to study how our results depend on the 
weights we assign to the sequence of the ten isochrones in 
Table[T]with the highest metallicities, which together repre- 
sent the 'thin disc'. As explained in Section [4.41 we consider 
three values for the observational scatter in [Fe/H]; 0, 0.1 
and 0.12 dex. 

GCS2 give the metallicity distribution of F and G 
dwarfs near the sun, which is a biased measure of the rela- 
tive numbers of stars formed with each metallicity because 
stars of a given mass but difi'erent metallicities have different 
lifetimes and luminosities and therefore probabilities of en- 
tering a magnitude-limited sample. Fortunately our model 
simulates this bias; we simply have to choose the weights 
of the isochrones such that the contribution of each chemi- 
cal composition to the modelled magnitude-limited sample 
agrees with the distribution of metallicities in GCS2. Let Wi 
be the proportion of the SFR which goes into stars of the 
ith chemical composition. Then the effective weight of this 
composition is Wi — FiWi/W , where W = FiWi and 

F.J,iB-V)J,r{^)^_^, (20) 

with (dA'^i/dr)dr the number of stars in the given colour 
and age range that we would have if the whole disc con- 
sisted of stars of the ith chemical composition. The Wi are 
chosen such that after convolution by an appropriate Gaus- 
sian distribution of measuring errors the effective weights Wi 
agree with the metallicity distribution determined by GCS2. 
The broken curve in Fig. 1111 shows the effective weights Wi 
that correspond to the intrinsic weights Wi of the isochrones, 
which are shown by the full curve - note that neither curve 
looks Gaussian because the bins have varying widths. We 
see that the effective weights are biased towards low metal- 
licities relative to the intrinsic weights Wi. 



0.25 ^ ^ , , , ^ , , , f , , , f , , , ^ , , r 




[Fe/H] 

Figure 11. The full curve joins the intrinsic weights Wi of the 
isochrones i = 3, . . . , 12 for the 'thin disc' component that yield 
the effective weights Wi of the isochrones in a magnitude-limited 
sample like that of GCS2 that are joined by the dashed curve. 
After convolution by a Gaussian distribution of observational er- 
rors with dispersion A[pg/jj] = 0.1 dex the Wi are consistent with 
the metallicity distribution of GCS2. 

The adopted metallicity distribution of the 'thin disc' 
component hardly affects the results of the best fit. All pa- 
rameters and the fit quality are very stable, as is shown in 
TableO As the construction of the low-metallicity compo- 
nent is not straight-forward, we also tested models with- 
out this component. The results are also presented in Ta- 
ble[51 For those models the age Tmax depends strongly on 
the metallicity distribution in the sense that the stronger 
the low-metallicity tail of the distribution, the larger is the 
recovered value of Tmax. Consequently, Tmax decreases as we 
increase our estimate of the errors AfcH in the GCS metal- 
licities because the larger the measurement errors, the nar- 
rower the distribution of intrinsic metallicities used in the 
model. The models without a low-Z component generally 
show better fit qualities than the two-component models, 
but broadly speaking all fits are in the same quality range. 
The old, low-metallicity component is constrained to a rela- 
tively small colour interval, whereas for the one-component 
model, the low-metallicity stars cover all ages and thus a 
wider colour interval, resulting in a better fit quality. All 
models with a significant low-metallicity component show 
very similar ages Tmax ~ 12.5Gyr. 

All the ages in Table [2] are larger than the value 
11.2 ± 0.75 Gyr obtained by BDBOO, principally because 
the reduction in the estimate of reddening shifts Parenago's 
discontinuity to the red, where the average age of stars is 
higher. The age could be brought down by lowering the he- 
lium abundance of the low-Z isochrones, which would shift 
the isochrones to redder colours. However, values of Y signif- 
icantly below the value of Yp from WMAP would be required 
to red uce Tmax appreciably. Interestingly, ICasagrande et aU 
(|2007l ) found very low Helium abundances for nearby K- 
dwarfs using Padova isochrones. 

For the following sections we decided to use two- 
component models with a 'thin-disc' metallicity distribution 
corresponding to measurement errors of Afch = 0.1 dex and 
an intrinsic dispersion of 0.14 dex. 
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Table 2. The parameters for the best fits to the data at different metallicity distributions and the fit quality x^. The metallicity 
distribution of the 'thin disc' component is characterised by the effective weights of the isochrones with the ten highest values of [Fc/H] 
in Table[T] 



^[Fo/h] 


effective weights Wi 


a 


13 


7 




n 


fio 




'thin disc' comp. 
















dex 


low high Z in % 






Gyr-i 


Gyr 


Gyr 


km 




with low met. tail 





(0.4, 3.5, 10.0, 16.0, 18.5, 19.5, 16.0, 11.0, 4.4, 0.7) 


-2.559 


0.356 


0.117 


12.557 


0.187 


55.187 


1.05 


0.1 


(0.0, 1.5, 8.0, 17.0, 22.0, 23.0, 17.0, 9.0, 2.5, 0.0) 


-2.549 


0.349 


0.117 


12.602 


0.149 


55.179 


1.04 


0.12 


(0.0, 0.7, 6.0, 17.0, 24.0, 26.0, 17.5, 7.5, 1.3, 0.0) 


-2.571 


0.350 


0.115 


12.601 


0.148 


55.232 


1.06 


without low met. tail 





(0.4, 3.5, 10.0, 16.0, 18.5, 19.5, 16.0, 11.0, 4.4, 0.7) 


-2.519 


0.385 


0.119 


12.606 


0.261 


57.157 


0.92 


0.1 


(0.0, 1.5, 8.0, 17.0, 22.0, 23.0, 17.0, 9.0, 2.5, 0.0) 


-2.673 


0.375 


0.121 


11.980 


0.201 


57.588 


0.94 


0.12 


(0.0, 0.7, 6.0, 17.0, 24.0, 26.0, 17.5, 7.5, 1.3, 0.0) 


-2.610 


0.376 


0.130 


11.782 


0.190 


57.975 


1.03 



Table 3. The effect of varying /3 on the fits. For 


a and Tmax 


upper 


limits of - 


-2.400 and 13.000 Gyr were set. 






a 




7 




Tl 


VlO 






fixed 


Gyr-l 


Gyr 


Gyr 


kms~^ 




-2.400 


0.250 


0.141 


13.000 


0.001 


48.257 


3.78 


-2.445 


0.300 


0.119 


13.000 


0.001 


52.615 


1.29 


-2.549 


0.349 


0.117 


12.602 


0.149 


55.179 


1.04 


-2.643 


0.420 


0.112 


12.281 


0.486 


57.449 


1.19 


-2.722 


0.500 


0.120 


11.605 


0.978 


58.963 


1.58 



5.2 Varying /3 

TableO shows the results of fixing l3, the exponent in the 
heating rate, and Fig. [12] shows the fits obtained with the 
best value (/3 — 0.349) and fits for values of /3 that are just 
too large (/? = 0.500) and clearly too smaU (/3 = 0.250) 
to be acceptable. For /3 = 0.250 the dependence of a on 
B ~ V is too flat and we obtain a poor fit to the number 
counts. Judging from the results, we are able to exclude 
values oi 13 < 0.28. The fit for /3 = 0.500 is relatively bad 
for the A'^-data of red stars and the corresponding a{B — 
V) is too flat for blue stars and too steep just before and 
around Parenago's discontinuity. In view also of the large 
X^, and the high velocity dispersion at birth (17.6 kms~^), 
we conclude that l3 ^ 0.50 can be excluded. 

When P is increased, the program increases vio and 
Tl in an effort to keep the general shape of o-{t) constant. 
Since the velocity dispersion of red stars strongly constrains 
vioT^ax! fhe increase in vio is compensated by a decrease in 
Tmax- As the isochroues have limted age ranges, we generally 
set an upper limit for Tmax at 13.0 Gyr. a and 7 are rela- 
tively stable for th e different values of / 3. a is only allowed 
to vary around the iKroupa et ahl l| 19931 ) value —2.7 within 
the interval (-2.4, -3.0). 

5.3 Varying the age 

We now examine the effect of optimising the other param- 
eters for flxed age Tmax. The favoured ages are higher than 
in BDBOO and surprisingly large in relation to the accepted 
age of the universe th = 13.7 Gyr. So we fixed Tmax at values 
^ 9 Gyr to see if acceptable younger models could be found. 



Table 4. The effect of varying Tmax on the fits. For the lowest 
ages considered, a was fixed at —2.400. 



a /9 7 Tmax Tl VlQ 

Gyr~^ fixed Gyr kms~^ 



-2.510 


0.338 


0.112 


13.000 


0.103 


54.758 


1.07 


-2.549 


0.349 


0.117 


12.602 


0.149 


55.179 


1.04 


-2.503 


0.357 


0.131 


12.000 


0.182 


55.495 


1.07 


-2.470 


0.363 


0.153 


11.500 


0.215 


55.499 


1.18 


-2.400 


0.375 


0.173 


11.000 


0.276 


55.696 


1.30 


-2.400 


0.383 


0.195 


10.500 


0.307 


55.951 


1.43 


-2.400 


0.393 


0.217 


10.000 


0.361 


56.124 


1.72 


-2.400 


0.404 


0.274 


9.000 


0.394 


56.594 


2.55 



Table|4] shows the results and Fig. [13] displays three fits to 
the data. 

7 increases with decreasing age as the ratio of old stars 
to young stars is an important constraint infiuencing both 
a{B — V) and dN/dC. a slightly increases with increasing 
7, but is relatively stable. As explained above, low ages are 
associated with higher values of (3, ti and vio- 

Judging from the values of listed in Table [4] ages 
above 11.5 Gyr are certainly accepted and ages of 10 Gyr 
and lower can be excluded. From Fig.[T3] we see that the 
problem with low ages lies in the reddest bins, in both the 
TV and the a data. The fit for 9 Gyr is definitely unable to 
represent the data and such low ages have to be excluded. 
For 10.5 Gyr, is still reasonable but the fit redwards of 
B — V = 0.55 looks bad. In view of the development of 
for ages between 11.5 and 10 Gyr, we can set a conservative 
lower age limit ~ 10.5 Gyr. 

5.4 The velocity dispersion at birth and after Tmax 

The correlations of the parameters ti and «io, which de- 
termine the velocity dispersion at birth and for the oldest 
stars, have already been discussed, vio only varies within a 
narrow range, as it is strongly constrained by the observed 
velocity dispersion redwards of Parenago's discontinuity. As 
Table[S] shows, the data allow a wide range of ti values that 
correspond to velocity dispersions at birth between and 
~ 15kms~^. The effect of this parameter is quite small as 
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55 




B-V 

Figure 12. Fits to the data for different values of /3. Upper panel: Total velocity dispersion, Lower panel: Projected number in the solar 
cylinder per colour bin. Green: Input Data, Long dashed, black: (3 = 0.250, Blue: /3 = 0.349, Short dashed, red: /3 = 0.500. 



Table 5. The effect of varying the parameter characterising the velocity dispersion at birth ri on the fits. 



O-birth O-fliial a /3 7 Tniax Tl Vw ^ 

kms~^ kms~^ Gyr^^ Gyr fixed kms~^ 



0.00 


58.90 


-2.400 


0.317 


0.117 


13.000 


0.000 


54.200 


1.11 


8.14 


59.46 


-2.519 


0.329 


0.116 


12.673 


0.030 


55.056 


1.09 


12.65 


61.82 


-2.549 


0.349 


0.117 


12.602 


0.149 


55.179 


1.04 


16.04 


61.36 


-2.635 


0.416 


0.119 


12.096 


0.500 


56.944 


1.17 


17.27 


62.22 


-2.728 


0.462 


0.119 


12.020 


0.800 


57.483 


1.39 



the power-law leads to a steep increase of velocity-dispersion 
within a time small compared to Tmax- As the dispersion is 
aheady 20kms^^ for the bluest bin, the adjustment of the 
parameters to fit a(B — V) is a minor issue. 

5.5 Studying the components of a 

As the results in Section[3]have shown, the three components 
of a show significantly different behaviours: their ratios vary 
with colour and Parenago's discontinuity is sharper in some 



components than in others. Therefore we now seek separate 
models for the evolution of each eigenvalue of . Since all 
three models have to satisfy the same A'^ data, and the errors 
on the eigenvalues of are larger than those on atotai , we fix 
the parameters that describe the star formation history at 
the best-fit values determined above for (Ttotai, namely a = 
—2.55, 7 = 0.117 and Tmax = 12.5 Gyr, and determined for 
each eigenvalue of only values of the heating parameters, 
/3, Tl and uio. The data for aw were represented by the 
polynomial fit shown in Fig.|4l 
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Table 6. Fits to the data for the three diflferent components of 
a. 





13 


n 


fio 








Gyr 


kms~^ 






0.307 


0.001 


41.899 


0.83 


ay 


0.430 


0.715 


28.823 


0.51 


aw 


0.445 


0.001 


23.831 


0.43 



Table[6] shows the best-fit values of the parameter and 
. The latter are quite small for ay and aw because the 
formal errors on these components are large. Compared to 
Aotai, Pu ~ 0.307 is lower, so the heating exponents for the 
other two components are higher than 0.349. Consequently, 
Pv = 0.430 and Pw = 0.445 are significantly larger than 
Pu- 

The three P values are consistent with Figure [5] in the 
sense that ai/a2 and cri/crs both decrease with increasing 
B — V and thus increasing mean age, and ai/cra decreases 
more steeply than (ti/(T2, indicating that Pw > Pv > Pu, 



Table 7. The effect of varying 7 at fixed a = —2.7 on the fits. 
As an upper limit for the age rmax wc apply 13.0 Gyr. 



a 


P 


7 


Tmax 


Tl 






fixed 




fixed 


Gyr 


Gyr 


km 




-2.700 


0.347 


0.140 


12.044 


0.147 


54.898 


1.21 


-2.700 


0.349 


0.120 


12.542 


0.144 


55.237 


1.07 


-2.700 


0.348 


0.100 


13.000 


0.136 


55.396 


1.06 


-2.700 


0.359 


0.080 


13.000 


0.163 


56.280 


1.12 


-2.700 


0.371 


0.060 


13.000 


0.194 


57.195 


1.47 



which is what we find. However, it has to be mentioned, that 
because of the rather large interval of acceptable values for 
P encountered in Section [5.21 we cannot exclude Pw ^ Pv- 

5.6 Different approaches to the SFR 

Q and 7 are strongly correlated. A high, positive value of 7, 
i.e. a higher SFR in the past, increases the number of red 
stars relative to blue stars. This effect can be cancelled by a 
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Table 8. The results for a SFR with two timescales. The two input parameters are displayed on the left, the output as usual on the 
right, followed by the characteristic model parameters i9 and p. Keep in mind the best for the chosen metallicity distribution was 1.04. 
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relatively flat IMF creating more blue stars. For the inves- 
tigation of different star form ation histories, we thus flxed 
a at the lKroupa et aO (Il993l ) value of —2.7 and tested dif- 
ferent values of 7. The best fit in Section [5.31 showed a less 
steep IMF, so we expect now to find that the model with 
lowest has a lower value of 7. TableQconfirms that this is 
the case. In all the above Sections, the SFR is decreasing: for 
the best fit from Section [5.31 it decreases by a factor 4.4 be- 
tween the beginning of star formation until now. Combining 
Tables [7] and m factors between 2.2 and 6.7 are plausible. 

For values of 7 ^ 0.1, we had to fix the age Tmax to 13.0 
Gyr for the reason given above. Our models generally need a 
large number of old stars. Since stars 10 Gyr and older move 
very slowly in the CMD, it is interesting to ask whether 
acceptable models with a reduced lower age limit can be 
obtained by permitting a high SFR early on. Moreover, a 
short early period of intense star formation is envisaged in 
the popular proposal that the thick disc for med as a result of 
a ma jor accretion event ~ 10 Gyr ago (e.g. IChiappini et al.l 
I200J ). 

We thus applied a SFR as described by equation (I12p . 
The timescale of the first formation epoch should be small 
compared to Tmax, so we tried A = 1, 2 and 3Gyr~^. The 
results for different parameters A are displayed in Table [8] 

The fraction p of solar-neighbourhood stars that be- 
long to the 'thick disc' is determined by A and A. Ta- 
ble|8j gives p for each model - observationally ljuric et al.l 
(|2008) found p ~ 0.11 for their definition of the thick disc, 
so we used this value as a guidance in adjusting the pa- 
rameters. We further characterise the scenarios by the ra- 
tio of the two terms at the beginning of star formation: 
1} = A exp[(A — 7)rmax]. Since models with very small A 
are indistinguishable from standard pure-exponential mod- 
els, we stopped lowering A when reached a value com- 
parable to that of the best-fitting standard model. Signifi- 
cantly, no two-exponential model achieved a lower x^ than 
the best standard model. 

For all scenarios, the parameters a, (3, t\ and wio are 
very stable; the disc-heating parameters do not depend sig- 
nificantly on the applied SFR, so there is no need for a fur- 
ther discussion. Compared to its values from Section [5.31 7 
has decreased for the lower ages to compensate for the in- 
fluence of the added 'thick disc' SFR term. The most strik- 
ing difference is that the favoured ages have decreased by 
~ 1.0 Gyr. The best results were achieved for A ^ 2.0 Gyr~^. 




Figure 14. The fits to the N-data for a standard SFR 10.5 Gyr 
model and a two-component (A = 2.0, A = 10~*) 10.7 Gyr model 
in the relevant B — V interval. 



For the SFR considered here, the lower age limit has to be 
lowered to ~ 10.0 Gyr. 

Fig. 1141 compares the fits to the A'' data for the 10.5 Gyr 
model of Section 15.31 and the 10.7 Gyr model with A — 
2Gyr~^ and A — 10~®. It shows how the intense first pe- 
riod of star formation improves the fits to the data of the 
red colo ur bins at lower disc ages. 

Roc ha-Pinto et al.l (2OOO) argue that the SFR of the 
solar neighbourhood is very irregular. They also find that 
the SFR shows an increasing tendency, which is incompat- 
ible with our findings. It is however interesting to ask if it 
is possible to achieve reasonable fits with a smooth SFR 
overlaid with factor vary ing in time according to Fig. 8 of 



overiaiQ witn lactor vary 
iRocha- Pinto et all (|2000l ) 



The best fit for this approach has 
a of 1.44 and is characterised by the following parameters 

a = -2.635 13 = 0.354 7 = 0.120 
Tmax = 12.266 Tl = 0.159 vio = 56.488, 

which are all in the range of acceptable values as determined 
above. The higher value of x^ results from additional fea- 
tures, which are produced by the varying SFR and which 
are incompatible with the data. We are, however, not able 
to exclude a SFR which is not smooth. As Figures [121 and [131 
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Figure 15. The effect of assigning a fiigti veiocity dispersion 
(Ttot = 85kms^^ to the low-metallicity component. The red line 
shows the best fit from above. 



show, there are features in the data that our models are not 
able to reproduce and which might be produced by epochs 
of enhanced star formation. 

We finally test the SFR described by eq uation (|13|) . 
which was proposed by lJust fc JahreissI (|2007l ). This SFR 
increases early on and then decreases. The best-fit model 
achieves = 1-15 with 



T2 « 14 Gyr 
Q = -2.400 
rn,ax = 12.686 



Tl 



0.334 
= 0.094 



b = 7.230 
vio = 54.392. 



In this model the SFR peaked ~ 10.3 Gyr ago and has since 
decreased by a factor ~ 4. The fit provided by this is of a 
similar quality as the best fits with the standard SFR, and 
the values of all comparable parameters are little changed 
from when we used the standard SFR. 



5.7 Different Approaches to the heating rate 

In our models we have so far assumed that the velocity dis- 
persion increases continuously with age. In the context of 
the Galactic thick disc, an obvious question to ask is what 



would be the effe ct of assigning a h igh velocity dispersion, 
ctot = 85kms~^ (jBensbv et al.ll2004l ). to the low-metallicity 
component. Fig. 1151 shows that this procedure produces an 
additional feature in the curve atot{B — V) that underlines 
again the narrow colour interval to which the low-metallicity 
component is confined. The reason for this is that the thick 
disk is not only metal-poor and old, bu t has a more complex 
structure jSchonrich fc Binnevl[2009bl ). 

In this context it is also interesting to study whether 
a saturation of disc heating, or a merger event producing a 
discontinuity in the time evolution of cr, or a combination 
of the two is compatible with the data. We adopt the al- 
ternative model (|17p of cr(r). This model introduces three 
additional parameters, Ti when heating saturates, T2 when 
there is a step increase in a and 77, the scale of that increase. 
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Figure 16. The fits for two different saturation scenarios with 
(^1, T2, r;) = (3Gyr, 9Gyr, 1.5) (red) and (5Gyr, 9Gyr, 1.5) 
(black, long dashed). The blue, short dashed curve is the best fit 
without saturation. The black arrows indicate the two additional 
discontinuities in the 3/9 Gyr scenario. 



Table[9] provides an overview of results obtained by fixing 
these parameters. 

We start by looking at a pure saturation of disc heating 
(7; — 1) with Ti = 5, 6 and 7 Gyr. Intuitively, this model is 
incompatible with the phenomenon of Parenago's disconti- 
nuity unless Ti > 7 Gyr, and that is what the results con- 
firm. For 77 = 1 we can firmly exc lude Ti ^ 5 Gyr. 

Consider now the proposal of lQuillen fc GarnettI (|200G| ) 
that Ti = 3 Gyr, T2 = 9 Gyr and tj > 1.5. Our best fit, 
obtained with rj = 1.75, is shown by the red curve in Fig. 1161 
The curve of a{B — V) moves from the high to the low side 
of the data around B — V — 0.55 as a result of two changes 
in slope at the points marked by arrows; at these points the 
main-sequence age corresponds to Ti and T2 . This figure and 
the high given in Table [§] exclude this model. 

To see whether the unwanted features associated with 
Ti and T2 could be washed out by an irregular star- 
formation history rather than the smooth one used in our 
models, we modelled the data using Ti = 3 Gyr, r2 = 9 Gyr, 
rj — 1.75 and overlaying our smooth SFR with a factor 
varyin g with time according to Fig. 8 of iRocha-Pinto et al.l 
(j2000h . The resulting fit to the data, shown by the long- 
dashed line in Fig. 1171 is even worse than that shown in 
Fig. 1161 Another possibility is that large errors in B — V 
(~ 0.023 mag) could smooth away the unwanted disconti- 
nuities. To test this hypothesis we folded a{B — V) for the 
(3 Gyr, 9 Gyr, 1.75) model with a Gaussian of 0.023 mag 
dispersion. Still the discontinuities did not disappear (short- 
dashed curve in Fig. ll7|l . 

Acceptable models can be found by allowing Ti to ap- 
proach T2 with corresponding adjustment to rj. For example, 
the fit provided by the (5 Gyr, 9 Gyr, 1.5) model is shown 
as the black curve in Fig. 1161 the fit to a{B — V) is clearly 
better than that for (3 Gyr, 9 Gyr, 1.75) model and nearly as 
good as the best fit without saturation (blue curve). Over- 
all one has the impression that models with Ti close to T2 
are merely approximating a power-law dependence with an 
appropriate step. 
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Table 9. The results for the assumption of a saturation of disc heating. The three input parameters are displayed on the left, the output 
as usual on the right. Keep in mind the best for the chosen metallicity distribution was 1.04. An upper age limit was set at 13.0 Gyr. 
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Figure 17. Attempts to smooth out the discontinuities in the 
(3 Gyr, 9 Gyr, 1.75) tr-plot (red curve) . A) Smooth SFR overlaid 
with a variation in time according to iRocha-Pinto et al.l (|2000| ) 
(black, long dashed). B) Folded with a Gaussian with a dispersion 
equivalent to the error in colour B — V (blue, short dashed) 



5.8 Why does a decrease redwards of the 
discontinuity? 

All plots of a(B — V), whether for models or data, show 
a counterintuitive decrease in sigma for the reddest bins: 
the velocity dispersion in our models increases with age, so 
decreasing a implies a decrease in mean age as one moves 
redwards of the discontinuity. Fig. [18] shows model age dis- 
tributions at B -V PS 0.50, 0.65 and 0.90. For the reddest 
colour we see only a smooth increase in numbers with age, 
reflecting the declining SFR, and a feature at high ages, 
resulting from the low-metallicty component. At the bluer 
colours we have strong features. From equation (|14p we see 
that they must reflect wider mass intervals yielding stars of 
the given colour at a certain time. The explanation for this is 
that in the vicinity of the turnoff, the isochrones run almost 
vertically and the colours of stars become nearly indepen- 
dent of mass for a significant range of masses. In Fig. [18] the 
blue age distribution for _B — 1/ « 0.65 shows eight peaks 
produced in this way, one from each of the eight metallicities 
used for the 'thin disk' component and also a low-metallicity 
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Figure 18. The age distribution of the 12.6 Gyr fit from Section 
l5.3l for three mean colours. 



feature at high ages. Because the dominant peaks all lie at 
ages higher than 5 Gyr, they raise the average age of stars 
in this colour bin above that of the bin for B — V — 0.91, 
which contains only stars too low in mass to have reached 
the turnofl. 



6 CONCLUSION 

We have updated the work of DB98 and BDBOO to take 
advantage of the re working of the Hipparcos catalogue by 
Ivan LeeuwenI |20o3), the availability of line-of-sight veloc- 
ities from the Geneva-Copenhagen survey, and significant 
improvements to the modelling. The latter include updated 
isochrones, better treatments of interstellar reddening, the 
selection function and the age-metallicity relation, and ex- 
ploration of a wider range of histories of star-formation and 
stellar acceleration. 

We have redetermined the solar motion with respect to 
the LSR. The new value (eq. 1) is very similar to that of 
DB98 but has smaller error bars. We have redetermined the 
structure of the velocity ellipsoid as a function of colour. 
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Again the results differ from those of DB98 mainly in the 
reduced error bars. 

The most striking change in results compared to BDBOO 
is an increase in the minimum age of the solar neighbour- 
hood Tmax • This iucrease is largely due to an improved treat- 
ment of interstellar reddening, which was overestimated by 
BDBOO. Models in which the SFR is oc exp(7T) fit the data 
best at Tmax ~ 12.5 Gyr and favour Tmax > 11.5 Gyr; these 
models yield a lower age limit of 10.5 Gyr. Models in which 
the SFR is a double exponential corresponding to the for- 
mation of the thick and thin discs favour ages in the range 
10.5 — 12 Gyr and yield a lower limit of 10.0 Gyr. For either 
model of the SFR, the lower age hmit of BDBOO, 9 Gyr, 
can be excluded. Our estim ates are in good agreem ent with 
the 12 Gyr age derived by ljust fc Jahreisd (|2003) and the 
11.7 Gyr obtained by SB09a. It is also in agreement with 
the individual stellar ages of GCS and GCS2, which include 
a significant fraction of ages between 10 and 15 Gyr. The 
study of Galacti c evolu tion in a cosmological context by 
[Hernandez et all ll200l[) y ields Tmax = 11 Gyr, while that 
of iNaab fc Ostrikeil ( 2006l l yields Tmax ~ 10 Gyr, our lower 
limit. Thus there is a considerable body of evidence that the 
solar neighbourhood started forming a remarkably long time 
ago. 

The strong correlation between the parameters a and 
7 t hat characterise the IM F and SFR, which was discussed 
bv iHavwood eUn (IT997l l and BDBO O, limit what w e can 
say about the IMF and SFR. We use a lKroupa et~all (|l993l ) 
IMF and find that the SFR is decreasing: the factor by which 
it decreases from the beginning of star formation until now 
is found to lie between ~ 6.5 and ^ 2.5. This co nclusion 
agrees with the findings of IChiappini et al.l (|l997f ) and of 
lJust fc Jahreisd (|2007), who used a non-exponential time 
dependence of the SFR, and what we find when the SFR is 
modelled by a sum of exponentials in time. 

Our conclusion regarding the time dependence of the 
SFR confiicts with the finding of BDBOO that the SFR was 
essentially fiat for a Salpeter IMF because the introduction 
of variable scale heights decreases the visibility of red stars 
relative to blue stars, so the models have to predict the ex- 
istence of a higher fraction of red stars than formerly. It 
also confiicts with the conclusion of SB09a that the SFR is 
only mildly decreasing for a Salpeter IMF because in their 
models many of the old stars in the solar neighbourhood are 
immigrants from smaller radii whi l e imm igration of young 
stars is negligible. [Hernandez et all (|200ll ) conclude that the 
SFR increases in the first ~ 4 Gyr and then d ecreases by a 
factor 2 until now, while iRocha-Pinto et al.1 (2000) argue 
that the Milky Way disc has a generally increasing but very 
irregular SFR. Neither picture is compatible with our re- 
sults. There ar e simi lar confiicts with the conclusions of 
iBertelli fc Na si' ('2001') that the SFR is broadly increasing 
and of lNaab fc Ostrikcr (2006) that the SFR increases early 
on and has been roughly constant over t he last ^ 5 Gyr. 

W e can exclude the scenario of lOuillen fc GarnetO 
l|2000l ') that disc heating saturates after 3 Gyr but at 9 Gyr 
the velocity dispersion abruptly increases by a factor of al- 
most 2. However, we are not able to exclude a later satura- 
tion of disc heating (after > 4 Gyr) that is combined with an 
abrupt in crease in dispersion more re cently than 9 Gyr ago. 
Similarlv. lSeabroke fc Gilmord l|2007l ) find that a saturation 
at > 4.5 Gyr is not excluded and that there is "extremely 



tentative" evidence of an abrupt feature in the age-velocity 
dispersion at 8 Gyr. Nevertheless, nothing in the data calls 
for early saturation of disc heating, and scenarios that in- 
clude it yield worse but formally acceptable fits to the data. 

Our favoured value /3 — 0.35 for the exponent that gov- 
erns the growth of utotai is in perfect agreement with the 
findings of GCS (0.34) and BDBOO (0.33). Moreover, the 
value of GCS3 (0.40) is also still in the range a llowed by our 
models. The classical value of ^ (|Wielenl|l977l ) yields rather 
bad fits and has to be regarded as the very upper limit for 
/3. 

For au we find /3 ~ 0.31, the same value as in GCS 
(0.31), but lower than that of GCS3, whose 0.39 is yet not 
out of range. For av we find (3 « 0.43, which is higher than 
the GCS and GCS3 values of 0.34 and 0.40. However, had 
they not ignored their oldest bins of stars, they would have 
obtained a larger value for /3. For aw we find /? « 0.45, 
in good agreement with GCS (0.47) and lower than the 
va lue of GCS3 (0.53), bu t higher than the 0.375 favoured 
bv lJust fc Jahreisd l|2007l ). 

The values and time dependencies of the ratios o"i/(T2 
and o"i/o"3 plotted in Fig. [5] provide important clues to the 
still controversial mechanism of stellar accelera tion. The 
original proposal (jSpitzer fc Schwarzschiidlll953D was that 
acceleration is a result of stars scattering off gas clouds. 
This process le ads to ch aracteristic axial ratios of the ve- 
locity ellipsoid. ISellwood (2008 ) has recently redetermined 
these ratios and finds ai/(T2 = 1.41 and ai/as — 1.61, both 
of which are smaller than Fig. [5] implies. This result sug- 
gests that scattering by spiral arms, wh ich increases a i and 
(72 but not (T3, plays an important role (|jenkinslll992t ). Fur- 
ther work is needed to quantify this statement in the light 
of Sellwood's recent work. 

Another indication that scattering by clouds is 
not alo ne responsible for acceleration is the finding of 
iHanninen fc Flvnnl ^2mi ) that in simulations of disc heat- 
ing by clouds, l3 = (0.21 ± 0.02), a value that we have 
excluded. For a disc heated by scattering off a combina- 
tion of gas clouds a nd massive black holes in the dark halo 
IHanninen fc Flvnnl (2002) find (3 = 0.42 and ai/as = 1.59, 
which is excluded by Fig. [5] 

It is widely believed that acceleration by spiral arms 
is responsible for the con centration of solar-neighbourhood 
stars in the (U, V) plan e (Raboud et al.' 1998; Dohnon '2000|; 
iDe Simone et al■ll2004^ ■ [Seabroke fc Gilmo re (2007) have ar- 
gued that this concentration undermines the concept of the 
velocity ellipsoid. However, it remains the case that the ve- 
locity dispersions au and ay increase systematically with 
age and it is not evident that any inconsistency arises from 
modelling this phenomenon as we do here just because spi- 
ral structure accelerates stars in groups rather than individ- 
ually. 

SB09a have recently argued that radial mixing plays an 
important role in disc heating and causes a higher value of 
P to be measured in the solar neighbourhood than charac- 
terises the underlying acceleration process. The origin of this 
effect is that stars migrate into the solar neighbourhood from 
small radii where the velocity dispersions are relatively high. 
Since the fraction of immigrants among the local population 
increases with stellar age, this effect enhances the velocity 
dispersion of old stars more than that of young stars, lead- 
ing to a larger effective value of /3. Our value for /3w is in 
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good agreement with the prediction of SB09a for a Hipparcos 
sample. 

The models of SB09a are more elaborate than ours not 
only in that they include radial migration but also in that 
they include chemical evolution. However, their parameters 
are determined by fitting to the GCS sample rather than the 
Hipparcos sample employed here. The current sample cov- 
ers a wider range of colours and is better defined than the 
GCS sample. Moreover, SB09a made no attempt to match 
the kinematics of the GCS sample but assumed the valid- 
ity of the description of the acceleration process given by 
BDBOO; they confined themselves to the metallicity distri- 
bution and Hess diagram of the GCS sample. Ideally one 
would fit simultaneously the local kinematics, metallicity 
and Hess diagram. However, given the impor t ance of radial 
migration implied by the studies of ' HavwoodI (|2008t ). SB09a 
and Schonrich & Binnov ( 2009b) , such a study would ideally 
include a more realistic treatment of the integrals of motion 
than the simple separability of the radial and vertical mo- 
tions assumed here and by SB09a. 

The cosm ic SF rate peaked at redshifts 1 — 2 (e.g. 
ILv et aLllioOTl ). which in the concordance cosmology corre- 
sponds to look-back times between 7.8 and 10.5 Gyr. These 
times are 1 — 2 Gyr later than the times Tmax at which star 
formation starts in our models, which is also the time at 
which the local star-formation rate peaked. It is also slightly 
later than the mean formation times, ~ 1.5±1 Gyr, of bulges 
in a recent series of larg e simulations of galaxy formation 
(jScannapieco et al.ll2009l ). On the other hand these authors 
found that the mean formation times of discs at _R ~ 10 kpc 
were 4±2 Gyr. The analogous time for the one of our models 
of the solar neighbourhood is 

_ _ _£-dtte^ 

' ""/.r'-dte.' 

The models listed in Table [2] yield values of r that range 
between 5.9 Gyr and 6.4 Gyr. T hese values fall at t he up - 
per end of the times obtained bv IScannapieco et al.l (120091 ). 
This overshoot may be connected to the fact that their mod- 
els have discs that are underweight by almost an order of 
magnitude; in the models disc formation may be artificially 
truncated. Thus our results are broadly in agreement with 
the results of ab-initio simulations of galaxy formation in 
the concordance cosmology, even though the age of the old- 
est solar-neighbourhood stars is remarkably large. At least 
some of these stars will be immigrants from small radii, and 
if they all are, star-formation will have started later than is 
implied by our values of Tmax. It is currently hard to place a 
limit on the fraction of the oldest stars that are immigrants 
because the models of SB09a assume that the disc's scale 
length does not increase over time, as is likely to be the 
case. 

Thanks to adaptive optics, gas-rich discs, in which stars 
are forming exceedingly rapidly, can now be studied obser- 
vationally at z ~ 2 (e.g. iGenzel et al.. 2008), soon after the 
oldest solar-neighbourhood stars formed. The discs observed 
are clumpy and highly turbulent, so it seems more likely 
that they will turn into bulges than a system as dynami- 
cally cold as the solar neighbourhood. Nonetheless, studies 



of these systems bring us tantalisingly close to the goal of 
tying together studies of 'galactic archaeology' such as ours 
with observations of galaxies forming in the remote past. 
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